from noise import *
from filters import *
import random

def e(l, i):
    return sim(i, l) + reg(l)
    
def sim(i, l):
    v_i = i.flatten()
    v_l = l.flatten()    
    
    r = v_i - v_l
    
    #print r[0 : 5]
    #print r.sum()
    #print multiply(r, r)[0: 5]
    #print multiply(r, r).sum()
    #print len(v_i)
    #print multiply(r, r).sum() / (len(i) * 1.)
    
    return multiply(r, r).sum() / len(r)
    
def reg(j):
    v_j = j.flatten()
    return absolute(gradient(v_j)).sum() / len(v_j)

# Imagenes en escala de grises
gray()    
clf()

img = nibabel.load('../../t1_acpc.nii.gz')
img = img.get_data()

z = 50

# Parte 1

#img_gaussian = gaussianFilter(img, 0.7, z)
#imshow(img_gaussian[:, :, z])
#show()

# Parte 2

#img_bilateral = bilateralFilter(img, 0.5, 0.3, z)
#imshow(img_bilateral[:, :, z])
#show()

# Parte 3

# Viendo cuanto da la energia
#print e(img, img)
#print e(img[:, :, z], gaussianFilter(img, 1.0, z)[:, :, z])

min = random.uniform(0.1, 0.5)
max = random.uniform(0.5, 1.0)

max_im = 30
max_it = 10

for i in xrange(max_im):
    sigma = random.uniform(0.1, 2.0)
    print sigma
    img_with_noise = gaussianNoise(img, min, max, sigma, z)   
    
    for j in xrange(max_it):
        sigma = random.uniform(0.1, 2.0)        
        img_gaussian = gaussianFilter(img_with_noise, sigma, z)
        print sigma, "\t", e(img_gaussian[:, :, z], img[:, :, z])
    
    print "\n"
    
    for j in xrange(max_it):
        sigma_1 = random.uniform(0.1, 2.0)
        #sigma_1 = 1.0
        sigma_2 = random.uniform(0.1, 2.0)
        #sigma_2 = 1.0
        img_bilateral = bilateralFilter(img_with_noise, sigma_1, sigma_2, z)
        print sigma_1, sigma_2, "\t", e(img_bilateral[:, :, z], img[:, :, z])
